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We demonstrate that models of resonant cyclotron radiation transfer in a strong field 
(i.e. cyclotron scattering) can account for spectral lines seen at two epochs, denoted SI 
and S2, in the Ginga data for GRB870303. SI, which extends 4 s, exhibits one line 
at ~ 20 keV, while S2, which extends 9 s, exhibits harmonically spaced lines at ~ 20 
and 40 keV. The midpoints of SI and S2 are separated by 22.5 s. Using a generalized 
version of the Monte Carlo code of Wang et al. (1988,1989b), we model line formation 
by injecting continuum photons into a static plane-parallel slab of electrons threaded by 
a strong neutron star magnetic field (~ 10^^ G) which may be oriented at an arbitrary 
angle relative to the slab normal. We examine two source geometries, which we denote 
"1-0" and "1-1," with the numbers representing the relative electron column densities 
above and below the continuum photon source plane. The 1-0 geometry may represent, 
e.g., a line- formation region levitating above the surface of the neutron star, or possibly 
a plasma-filled flux tube illuminated from below. The 1-1 geometry, on the other hand, 
corresponds to line formation in a semi-infinite atmosphere at the surface of a neutron 
star. We apply rigorous statistical inference to compare azimuthally symmetric models, 
i.e. models in which the magnetic field is parallel to the slab normal, with models having 
more general magnetic field orientations. If the bursting source has a simple dipole field, 
these two model classes represent line formation at the magnetic pole, or elsewhere on 
the stellar surface. We find that the data of SI and S2, considered individually, are 
consistent with both geometries, and with all magnetic field orientations, with the 
exception that the SI data clearly favors line formation away from a polar cap in the 
1-1 geometry, with the best-fit model placing the line- forming region at the magnetic 
equator. Within both geometries, fits to the combined (S1-I-S2) data marginally favor 
models which feature equatorial line formation, and in which the observer's orientation 
with respect to the slab changes between the two epochs. We interpret this change as 
being due to neutron star rotation, and we place limits on the rotation period. 

Subject headings: radiative transfer — line: formation — gamma rays: bursts 



1. Introduction 

The low-energy spectral lines detected between ~ 20-50 keV in the spectra of some gamma-ray 
bursts (GRBs) provide the most powerful means by which an analyst can probe the burst source 
environment. This is because the study of lines, unlike the study of continuum spectra, can yield 
precise information about the temperature and column depth of the line-formation region, the 
velocity and orientation of that region, and the strength and orientation of any magnetic field in 
that region. 

Candidate lines in GRB spectra were first reported by Mazets et al. (1980, 1981), who detected 
single absorption-like dips and troughs in the spectra of ~ 15% of bursts observed by the Konus 



detectors of Venera 11 and 12. Hueter (1987) reported similar absorption-like features in the spectra 
of two bursts observed by the HEAO-1 A4 detector, and the Ginga Gamma- Ray Burst Detector 
(GBD) observed three bursts-GRB870303, GRB880205, and GRB890929-whose spectra exhibited 
harmonically spaced absorption-like line candidates with centroid energies between ~ 20-50 keV 
(Murakami et al. 1988; Fenimore et al. 1988; Yoshida et al. 1991; Graziani et al. 1992, 1993; Freeman 
et al. 1999a, hereafter Paper I). A spectrum from an earlier epoch of GRB870303, denoted SI, was 
subsequently found to exhibit a single absorption-like line candidate at ~ 20 keV (Graziani et 
al. 1992, who used S2 to denote the other line candidate spectrum of GRB870303; Graziani et 
al. 1993, Paper I). Fits with phenomenological continuum-plus-line(s) models demonstrate that the 
Ginga line candidates have statistical significance ~ 10~^-10~^. The enhanced low-energy coverage 
of the Ginga GBD (E'low ~ 1-5 keV, compared with E'low k, 20 keV for Konus and HEAO-1 A4), 
allowed analysts to rule out the possibility that these line candidates could be explained by a sudden 
change in continuum shape. 

Mazets et al. posited that the Konus line candidates were formed by either cyclotron absorption 
or cyclotron scattering in the strong magnetic fields of galactic neutron stars (~ 10^^ G). The 
observation of harmonically spaced line candidates gave further credence to this hypothesis, because 
the quantization of an electron's energy perpendicular to a strong field can lead to the formation 
of evenly spaced lines with energies En ~ 11.6ni?i2 keV. Rigorous proof of the viability of this 
hypothesis came when a number of analysts used theoretical models of cyclotron scattering (and not 
cyclotron absorption) to produce emergent spectra from line-forming regions, which were found to 
compare favorably with the Ginga data (Wang et al. 1989a; Alexander & Meszaros 1989; Nishimura 
et al. 1992; Wang, Wasserman, & Lamb 1993). For instance, Wang et al. (1989a) used a Monte 
Carlo radiative transfer code (Wang, Wasserman, & Salpeter 1988, 1989b; Lamb et al. 1989) to 
demonstrate that the harmonic line candidates in the GRB880205 spectrum could be produced in a 
line-forming layer with magnetic field 1.7 x 10^^ G and column density 1.2 x 10^^ cm^^. Freeman 
et al. (1992) later demonstrated that deepening the scattering region behind the continuum source 
led to a substantially better fit to these same data. 

The cyclotron scattering interpretation for line candidates went largely unchallenged when it 
was consistent with the prevailing view that GRB sources were neutron stars residing in a thick 
disk in the Milky Way (with scale height ~ 2 kpc; see, e.g., Higdon & Lingenfelter 1990, Harding 
1991). However, recent developments, while not challenging the cyclotron scattering picture per 
se, have led many to question whether the reported line candidates actually exist at all. First, the 
discovery of optical transients (OTs) associated with GRBs (e.g. van Paradijs et al. 1997 and refer- 
ences therein), and the apparent determination of redshifts for five of them— GRB970508 (Metzger 
et al. 1997), GRB971214 (Kulkarni et al. 1998), GRB980613 (Djorgovski et al. 1999), GRB980703 
(Djorgovski et al. 1998), and GRB990123 (Kelson et al. 1999)— have indicated that some (if not all) 
GRBs occur at cosmological distances. Second, there have been no reports of definitive detections 
of line candidates in spectra of the Burst and Transient Source Experiment Spectroscopy Detectors 
(BATSE SDs) on the Gompton Gamma-Ray Observatory (Palmer et al. 1994, Band et al. 1996, 



Briggs et al. 1996, 1998). These developments, combined with the apparent difficulty of forming 
lines in the circumstellar environments of cosmological burst sources (cf. Stanek, Paczyhski, & 
Goodman 1993 and Ulmer & Goodman 1995, who attempt to account for lines by invoking gravi- 
tational femtolensing) , have led many to adopt the viewpoint that all bursts are cosmological and 
that the reported line candidates, for whatever reason, are not real. 

This viewpoint may be intuitively reasonable on its surface. However, while we do not seek 
to disprove the posited cosmological origin of the GRBs listed above, we feel that there are several 
reasons why we should continue to test the cyclotron scattering model. While none are completely 
compelling by themselves, they cumulatively indicate that the complete solution of the GRB mys- 
tery may have not yet been provided by the discovery of GRB redshifts. 

The first reasons are the suggestive pieces of direct observational evidence that there are 
two or more classes of classical gamma-ray bursts. These pieces of evidence range from the well- 
established to the speculative. Well-established is the identification of two classes with differing 
light-curve morphologies (Kouveliotou et al. 1993, using burst duration; Lamb, Graziani, & Smith 
1993, using light curve variability). Further studies lend credence to the reality of separate classes 
(Kouveliotou et al. 1996; Katz & Canel 1996), but it is not yet proved that the correlation of burst 
properties cannot be explained using a single underlying mechanism. 

Less well-established, but still highly suggestive, evidence is provided by burst repetition 
(e.g. Quashnock &: Lamb 1993, Wang & Lingenfelter 1995). Unlike classes based on light-curve 
morphology, repetition indicates directly that there must be a galactic component to the overall 
GRB source population, since it is generally considered impossible for cosmological bursters caused 
by one-time events such as neutron star-neutron star or neutron star-black hole mergers to repeat. 
Evidence for repetition is absent from the Third BATSE (3B) catalog (Meegan et al. 1996, Tegmark 
et al. 1996), but Graziani & Lamb (1996) have questioned how the systematic errors of this catalog 
are computed. (For details on the BATSE burst location algorithm, see Pendleton et al. 1999). 
Also, in October 1996, after the publication of the 3B catalog, BATSE observed four bursts in less 
than two days which all came from directions consistent with a single source. Graziani, Lamb, &: 
Quashnock (1998) perform simulations which indicate that the probability of such a spatial and 
temporal coincidence of four distinct bursts from four different sources is 3.1 x 10^^. The proba- 
bility is increased to 1.6 x 10^^ if only three bursts occurred, but one of the three would have to 
be the longest ever observed by BATSE. 

The most speculative piece of evidence is the apparent lack of consistency that recently localized 
bursts show when examined in other wavebands. For instance, since the beginning of 1997, there 
have been eleven GRBs observed with the BeppoS AX WFC and/or the RXTE ASM detectors whose 
positions have been refined through the use of the Interplanetary Network (IPN).^ Of these eleven 
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bursts, X-ray afterglows, OTs, and radio afterglows were detected for eight, six, and three of them 
respectively. (Two have measured redshifts and are thus conclusively cosmological— GRB971214 
and GRB990123.) The lack of counterparts may be, of course, more simply explained by invoking, 
e.g., selection effects, rather than by invoking a separate class of bursters whose broad-band behavior 
is such that they do not show optically transient emission. Also, if there is a separate class of 
bursters which are defined by counterpart behavior, it would not necessarily have to be associated 
with the Milky Way. 

Another reason for us to continue to test the cyclotron scattering model is the result of the 
exacting statistical analysis of the 3B catalog performed by Loredo & Wasserman (1995, 1998a, b). 
They determine that the 3B catalog is consistent with the hypothesis that there is a component of 
the GRB source population residing within the galaxy, with that component comprised of either 
dim local halo sources at distances ^ 1 kpc (accounting for up to ^ 60% of all bursts), or luminous 
halo sources at distances ^ 50 kpc (accounting for up to f» 10% of all bursts). Studies show that 
neutron stars which receive high initial kick velocities when they are formed can populate the halo 
(e.g. Li & Dermer 1992; Podziadlowski, Rees, & Ruderman 1995; Bulik, Lamb, & Coppi 1998). Such 
"hybrid" models with both galactic and cosmological bursters fit to the SB data better than models 
with only cosmological bursters, in part because the existence of a highly isotropic, extragalactic 
burst component weakens greatly the isotropy constraints on any anisotropic Galactic component. 
However, the fit is not so much better as to rule out models containing only cosmological bursters. 

A final reason to continue to test the cyclotron scattering model is the fact that the BATSE 
SDs may lack the low-energy spectral sensitivity that is necessary for low-energy lines to be easily 
detected. The gain settings on the individual SDs differ; those with the highest gain settings can, in 
principle, observe GRBs at energies ^ 10 keV. An electronic artifact discovered after launch affects 
energy calibration such that spectra are distorted in the first ~ 10 channels above the low-energy 
cutoff (the so-called "SLED" effect; see Band et al. 1992). Despite this, studies using simulated 
Ginga line candidate spectra indicated that the SDs were still capable of detecting low-energy 
line candidates (Band et al. 1995), and no line candidates were definitively detected during initial 
visual searches of those spectra with the largest signal-to-noise ratios (Palmer et al. 1994, Band 
et al. 1996). An automated line candidate search algorithm designed by the BATSE SD team 
(Briggs et al. 1996) was then applied to spectra in 117 bright bursts for which there is at least one 
spectrum with signal-to-noise > 5 at ~ 40 keV (Briggs et al. 1998). This automated search, which 
is considerably more sensitive than a visual search, yielded 12 candidate spectral line candidates for 
which the change in x^ between the continuum and continuum-plus-line fits is > 20 (significance 
< 5 X 10^^). Perhaps problematically, all the line candidates except one are emission-like lines 
observed at ~ 40 keV, with the one exception being an absorption-like line candidate observed at 
~ 60 keV. Briggs et al. (1998) cannot rule out the possibility that sharp breaks in the continua, 
instead of lines, cause the observed features, but such breaks would be inconsistent with low-energy 
Ginga GRB data. Briggs et al. estimate the ensemble chance probability of the most-significant 
feature as ^ 10"'^, and state that few of these features, if any, result from statistical fiuctuations. 



However, these lines should not be considered definitively detected, as the contemporaneous data 
from other SDs is still being examined (Briggs et al. 1999). 

We feel that the aforementioned reasons give us ample justification to test the cyclotron scat- 
tering model further. In this paper, we use it to examine the lines^*^ exhibited by GRB870303 
SI and S2. Successful fits of the cyclotron scattering model to these data, when paired with the 
successful fits of cyclotron scattering models to the data of GRB880205 by Wang et al. (1989a) and 
Freeman et al. (1992), would further strengthen the evidence supporting our contention that that 
some (though not all) GRBs are galactic in origin. 

In §2.1, we describe the spatial geometry of the line-formation region. We assume it to be a 
static plane-parallel slab of electrons threaded by a strong magnetic field, which can be oriented 
at an arbitrary angle relative to the slab normal (unlike in Wang et al. 1989a, where the field was 
parallel to the slab normal). The assumption of a static slab is strictly valid only for bursters located 
within ~ 100 pc, if line formation occurs at the magnetic pole, since otherwise the inferred burst 
luminosity would exceed the Eddington limit and radiation forces would blow the line-forming layer 
off in a wind (see Isenberg, Lamb, & Wang 1998a, who modify the code to examine the theoretical 
ramifications of a wind). As noted above, such distances are consistent with current observational 
limits (Loredo & Wasserman 1998b). If line formation occurs away from the magnetic pole, then 
the assumption of a static slab may be valid even for luminous halo bursters at distances ~ 100 kpc 
because of the confinement provided by closed magnetic field lines (e.g. Zheleznyakov and Serber 
1994, 1995). We examine two slab geometries, which we denote "1-0" and the "1-1", where the 
numbers represent relative electron column densities above and below the continuum photon source 
plane. A slab illuminated from below represents a line-forming region in the magnetosphere of a 
neutron star, while a source plane embedded within a slab corresponds to line formation within a 
semi-infinite neutron star atmosphere. 

In §2.2, we summarize the physics of radiation transfer in strong magnetic fields that is incor- 
porated into the Monte Carlo code we use to generate spectra (Wang et al. 1988, 1989b; Lamb et 
al. 1989; Lamb, Wang, & Wasserman 1990; and references therein), and in §2.3 we provide exam- 
ples of these spectra. As noted above, resonant cyclotron scattering, and not cyclotron absorption, 
describes the peculiar properties of the Ginga lines— the comparable strengths of first and second 
harmonics, the absence of third and higher harmonics, and the narrowness of the lines. The appear- 
ance of the lines is greatly affected by Raman scattering, i.e. resonant cyclotron scattering in which 
the electron is excited from the ground state directly to the second or higher Landau level, but then 
deexcites back to the ground state indirectly via intermediate energy levels. In Raman scattering, 
the original photon is destroyed, and two or more photons are created (or spawned) as the electron 
deexcites. Spawned photons with energies near that of the first harmonic line alter its profile and 



^"The rigorous statistical analysis carried out in Paper I on the data of GRB870303 SI and S2 demonstrates that 
when the data are considered jointly, the significance of the continuum-plus-lines model is ~ 10^*. Thus we feel that 
we may drop the word "candidates" when referring to these lines throughout the remainder of this paper. 



reduce its equivalent width to a value similar to that of the second harmonic line. Because the 
majority of photons which undergo scattering at the second and higher harmonics are destroyed, 
the lines have an approximately absorption-like profile. The third harmonic is not seen because 
its optical depth is small compared to that of the second harmonic (rs ~ 0.05r2). Narrow lines 
result from the fact that the line-forming region is optically thin to all continuum photons except 
for those with energies equal to the first and second harmonic energies. Thus scattered photons 
at ~ 20 keV, and not continuum photons ^ 1 MeV, dictate the temperature of the line-forming 
region, and the line width. 

In §3 we summarize the statistical concepts which we use in this paper. These concepts are 
described in more detail in Paper I and Freeman et al. (1998b; hereafter Paper III). In those works, 
we present general, rigorous methodologies that address the problem of establishing the existence of 
a line in a spectrum, that are based upon both the so-called "frequentist," and Bayesian, paradigms 
of statistical inference. In this work, instead of establishing the existence of lines, our statistical goal 
is to compare azimuthally symmetric models of line formation, i.e. models in which the magnetic 
field is parallel to the slab normal, with models having arbitrary magnetic field orientations. If the 
GRB source is a neutron star with a simple dipole field, these two classes may be interpreted as 
representing line formation at the magnetic pole, or elsewhere. We use a rigorous method of model 
comparison (described in §3.2) to select between polar cap and more general models, demanding, 
e.g., that the increase in the quality of the fit of general models be sufficiently great to justify 
considering of the ramifications of line formation away from the pole. Once best-fit models are 
selected, we use a rigorous method of parameter estimation (described in §3.3) to place limits on 
each of the model parameters. 

We apply the cyclotron scattering model to the data of GRB870303 SI and S2 in §4. SI, 
showing 4 s of data, exhibits a saturated line at ~ 20 keV, while S2, showing 9 s of data, exhibits 
two harmonically spaced lines at ~ 20 and 40 keV. The midpoints of the SI and S2 time intervals 
are 22.5 s apart. In Paper I, we establish the evidence for the GRB870303 lines using simple 
phenomenological models. In that work, we establish the frequentist significance of, and Bayesian 
odds favoring, the SI line to be 3.6 x 10^^ and 114:1, respectively; for S2, the respective figures 
are 1.7 x 10"'' and 7:1, while for the combined (S1-FS2) data, they are 4.2 x 10"^ and 40,300:1. In 
joint fits to the combined (S1-I-S2) data in this work, we find that the best-fit values of fi and/or (j) 
change during the 22.5 s between SI and S2 for both the 1-0 and 1-1 geometries, meaning that the 
orientation of the observer changes with time. We interpret the change in orientation by invoking 
neutron star rotation, and in §4.5 (and the Appendix) we describe how we place limits on the 
rotation period. 

In §5, we discuss our results. 



2. Cyclotron Scattering in Strong Magnetic Fields 

In this paper, we assess the hypothesis that the hnes exhibited at ~ 20 keV in GRB870303 
SI and ~ 20 and 40 keV in GRB870303 S2 were formed within a strong magnetic field (B ~ 
10^^ G). We make the further assumption that the lines were formed on the surface of, or near, a 
neutron star, the only astronomical object where such strong field strengths have been observed. 
We assess the hypothesis by generating spectra with a Monte Carlo code that numerically treats 
radiation transfer in strong fields (Wang et al. 1988), and fitting these spectra to the observed data. 
In this section, we describe the spatial geometry of the line-forming region, and then summarize 
the physics of radiation transfer included in the Monte Carlo code. We then provide examples of 
spectra generated with the code, so as to build the reader's intuition before we discuss the results 
of model fitting in §4. 



2.1. Spatial Geometry of the Line-Forming Region 

We assume that line formation occurs in a plane-parallel slab of electrons, which has infinite 
horizontal extent and height similar to a neutron star atmospheric scale height {h ~10~^i?Ar5). 
This slab contains an electron-proton plasma, threaded with a magnetic field B which is oriented 
at angle \I' with respect to the slab normal h (Figure 1).^^ We neglect variations in the magnitude 
and direction of B within the slab. 

Continuum photons are injected into the slab at a source plane, travel through it, and emerge 
from one of its faces. We refer to photons that emerge from the top of the slab as "transmitted" and 
those that emerge from the bottom as "refiected." Transmitted photons reach the observer, whose 
orientation relative to the slab is given by 6, the polar angle from the slab normal (or equivalently, 
/i = cos^), and (j), the azimuthal angle from the projection of B onto the slab (Figure 1). 

The location of the source plane with respect to the slab determines the geometry of the 
system. We apply the nomenclature "1-x" to denote geometries; 1 and x (a number) represent the 
relative values, above and below the source plane, of the electron column density N^ between the 
photon source plane and the observer. We examine two geometries in this work (Figure 2). The 
1-0 geometry is a slab illuminated from below; this geometry is used by Wang et al. (1989a). This 
geometry represents a line-formation region physically separated from an isotropically-emitting 
source of continuum photons, e.g. an illuminated flux tube. We assume that the reflected photons 
return to the neutron star surface, where they are thermalized, i.e. absorbed by non-resonant inverse 
magnetic bremsstrahlung. We use a 1-1 geometry to model line formation in an isothermal semi- 
inflnite neutron star atmosphere. Slater et al. (1982) and Wang et al. (1988) determine that the 
mean number of scattering events between the source and the top edge that a resonant photon 



^Wang et al. (1989a) consider only * = in their fits to the Gmga data of GRB880205. 



experiences prior to its escape approaches a limiting value as the atmosphere becomes semi-infinite 
(i.e. a l-oo geometry). The number of scattering events experienced in the 1-1 geometry is within 
~ 10% of its limiting value, while the number experienced in the 1-4 geometry (used, e.g., in 
Freeman et al. 1992) matches the limiting value (Isenberg, Lamb, & Wang 1998b). We do not 
use the 1-4 geometry, despite its greater accuracy, because the overall number of scattering events 
(both resonant and non-resonant) is approximately three times that using the 1-1 geometry; the 
computer time required for the 1-4 simulation increases proportionally. Also, in the 1-1 geometry, 
the reflection symmetry of the line- forming region allows us to limit injection of photons to a 
hemisphere facing the observer. The sum of the resulting transmitted and reflected spectra is 
equivalent to the spectrum which would result from spherically isotropic photon injection. By 
using hemispherical input, we reduce the amount of computer time per simulation by a factor of 
two. 



2.2. The Physics of Resonant Cyclotron Scattering 

Here, we summarize the physics of strong field radiation transfer that is incorporated into the 
Monte Carlo code. Unless we specify otherwise, the reader may find more details on any of the 
elements of the code described below in Wang et al. (1988). 

We assume that the line-forming plasma is cold: A;Tg y ^ Eb = ^ B k, ;^^i3, where 7 is the 
Lorentz factor (~ 1 for the electron energies assumed in this work), and the symbol || indicates that 
it is only the energies parallel to B that follow the (continuous) Maxwell-Boltzmann distribution. 
The allowed electron energies perpendicular to the field are the Landau levels. The spacing of these 
levels is assumed to be much larger than typical electron energies, so that they are not collisionally 
populated, and the photon densities are assumed to be small, so that they are also not radiatively 
populated. Thus, in each scattering event, the initial and final electron energy state is the Landau 
ground state (n = 0). 

Electrons which are excited to states n > 2 may either deexcite directly to the ground state 
(in which case the photon has undergone resonant scattering), or it may reach the ground state 
via intermediate excitation levels [Raman scattering). Raman scattering is dominant in weak fields 
(■g- <C 1, where Be is the critical field strength, 4.4x10^^ G, at which the first harmonic energy 
matches the electron mass energy 511 keV). For instance, for the magnetic field strengths considered 

in this paper, electrons excited to n = 2 have probability ~ 1 % = 1— -g- ~ 0.95 of deexciting 

to the first excited state (Daugherty & Ventura 1977). When this occurs, the incident photon with 
energy E ~ 2Eb is destroyed, and two photons with energies E ~ Eb are spawned. The second 
and higher harmonics thus have an absorption-like line profile in weak fields, while the spawned 
photons act to reduce the equivalent width of the first harmonic line from the value it would have 
had if Raman scattering did not occur. Because the equivalent widths of the second and higher 
harmonic lines are oc (-^)"~^, the code does not treat excitation to states n > 4. 
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The code uses scattering cross-sections provided by Herold (1979) and Daugherty & Ventura 

(1977). Herold states how to formulate exact relativistic scattering cross-sections, and provides 

leading-order expressions for the case where the initial and final electron states are the Landau 

ground state, valid in the limit 

E V Be 

<1. 1 



The code uses these expressions, which are valid in the energy and magnetic field strengths regimes 
which are of interest in this work. However, these expressions make no provision for the natural line 
width of the resonances, and as a consequence the scattering amplitudes diverge at the resonances. 
To correct the divergence, the natural line- width, r„ = n-^ — ^, is added (Wasserman & Salpeter 
1980, Harding & Daugherty 1991, Graziani 1993; a is the fine-structure constant, and the line-width 
is given in units of energy). The code also uses second and third harmonic resonant scattering cross- 
sections derived in the non-relativistic limit by Daugherty & Ventura, who use Fermi's Golden Rule. 
The first harmonic cross-section given by Daugherty & Ventura does not as accurately portray the 
line wing profile, hence the Herold cross-section, with finite line-width added, is used instead. 

A problem with this "hybrid" procedure which mixes exact relativistic resonant cross-sections 
and non-relativistic higher harmonic cross-sections is that in the treatment of electron-photon 
scattering, there is the implicit approximation 



E' 



oo^o + oo^i-^o|^ + yi \ai\^ , (2) 






where a^ is the matrix element for the i^^ scattering channel. This approximation is exact for 
scattering at the first harmonic both in the line core (within a few Doppler widths of line center) 
and wings (far from the line center). For higher harmonic scattering, this approximation is only 
valid near the line center. If there is significant scattering at energies between the first and second 
harmonics, then the matrix cross-term which is missing from this expression may be important. 
Consequently, we limit use of the code to modeling slabs which are not optically thick at the first 
harmonic line wings (see below; also Wasserman &: Salpeter; Lamb et al. 1989): 

where ai is the ratio of the natural line width to the Doppler energy width for the first harmonic, 

1 2kT I 

Wth is the thermal electron velocity along the magnetic field lines (= a/ ^° ), Ng is the column 
density of electrons in the slab, and ri is the polarization-, angle-, and frequency-averaged optical 
depth in the first harmonic. For relevant temperatures (kT^u '^ 5 keV, see below), Ne ^ 10^^ cm~^. 

(If Ne ^ 10^^ cm~^, the line- forming slab may be moderately optically thick to photons at the 
harmonic energies but will be optically thin to continuum photons. Thus scattered photons at ~ 
20 keV, and not continuum photons ^ 1 MeV, dictate the temperature of the line-forming region. 
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Lamb et al. 1990 derive what they dub the Compton temperature, Tq [which equals Tg n in this 
work, but not generahy], for this density regime, and find that [1] slabs reach this temperature on 
timescales t ~ 10~^ s, and [2] kTc is a fraction [~ 25%] of the first harmonic energy. Thus the hues 
are narrow. This result demonstrates the appropriateness of the originally stated assumption that 
lines are formed in a cold plasma.) 

The resonant cross-sections are averaged over the initial polarization states and summed over 
the final states. Polarization-averaged cross sections are appropriate for first harmonic scattering in 
optically thick media when the vacuum contribution to the dielectric tensor dominates the plasma 
contribution, a condition we may state as ( iq^"" -3 ) (tot^tt) ^ 1, where He is the electron number 
density. This condition holds for the physical conditions which we explore here, and is, in fact, less 
limiting than the requirement that line- forming regions be optically thin in the line wings. 

The line profile at each harmonic is created by making a Lorentz transformation to the lab 
frame, and averaging over /(p), the one-dimensional electron momentum distribution (along B). 
The total scattering profile is assumed to be the sum of the profiles for the individual Landau levels 
(including the continuum contribution). Profiles are averaged over azimuthal angle, (/>, despite the 
fact that for ^ 7^ 0, azimuthal symmetry is broken. The (?i>-dependent part of the scattering cross 
section is only significant in the line wings and continuum, and in the present work we considerly 
only line-forming regions which are not optically thick in the line wings. 

While we have used the terminology "line core" and "line wings," whose meanings are intu- 
itively well known, we can define them mathematically as fulfilling the conditions |^^| ^ 1 (core) 
and \xn\ 3> 1 (wings), where hb is the cosine of the angle between B and the direction of propaga- 
tion of the photon, and Xn = 7"th\ is a dimensionless energy shift (with the denominator being 

the Doppler width associated with the n harmonic energy En = tlEb', Wasserman & Salpeter). In 

the line core, the thermal electron distribution dominates the profile so that it is oc exp 

while in the wings, the tail of the Lorentzian distribution dominates so that the profile is oc a„a:^^, 
where a„ is again the ratio of the natural line width to the Doppler energy width. If ^b = 
the line wings extend to x{^b) = 0, i.e. the line profile is purely Lorentzian for all x. This is 
because the Doppler effect vanishes, to first order in ^, for photons propagating perpendicularly to 
B. Wasserman &: Salpeter showed that for the first harmonic, the core-wing boundary appears at 
-^ ~ 2.62 — 0.191og(-r^). We refer to the wings at energies below and above the line center as 

Ms °^lMsr ° 

the red and blue wings respectively. 

Lamb et al. (1989) showed that relativistic kinematics has a significant effect on the shape of 
the absorption profile, even in the limits E, kT^n <^ nieC^. For zero natural line width, relativistic 

kinematics prohibits scattering at the n harmonic above a cutoff energy Ec = " ^ ''^ -, 

where 6 = -^ (Daugherty & Ventura 1978; Harding & Daugherty 1991; see the Appendix of Wang 
et al. 1993 for a physical derivation). Wasserman &: Salpeter show that for the physical conditions 
where electron recoil is important, photons will escape more readily in the red wing than in the 
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blue. This result motivated Lamb et al. (1989), Wang et al. (1989b), and Wang et al. (1993) to 
make the simplistic assumption that the resonant scattering profile was zero for E > Ec, i.e. to 
ignore the effects of the finite natural line width. For E < Ec, the effects of finite line width were 
included. The simplification led to the appearance of spikes in some cyclotron scattering spectra at 
energies just above that of the first harmonic (e.g. as ^ ^ for ^ = 0). This spike contained both 
photons which were scattered at the first harmonic, and spawned photons resulting from Raman 
scattering at the second and higher harmonics, which then immediately escaped from the slab. The 
code used in this work does include the effect of finite natural line width for E > Ec so that the 
resonant scattering profile is now small but finite above Ec- With this enhancement, the spikes are 
smeared by scattering and no longer appear. We stress that even though the scattering profile is 
finite above the cutoff energy when the effects of natural line width are properly treated, the profile 
still falls off sharply above Ec, leading to a strongly asymmetric line shape in some spectra. 



2.3. Examples of Monte Carlo Spectra 

To help build the reader's intuition before describing the results of spectral fits to the data of 
GRB870303 SI and S2, we present examples of spectra produced with the Monte Carlo code. In 
Figures 3 and 4, we show output counts spectra produced with ^ = 0, for the 1-0 and 1-1 geometries 
respectively. In Figure 5 we show spectra produced with ^ = ^, for the 1-0 geometry (the 1-1 
spectra are similar and are not shown). Each spectrum is produced assuming B12 = 1.7 and A^e,2i 
= 0.6. The relative strengths of the harmonics in each figure may be understood as arising from 
the interplay of line-of-sight column density (i.e. whether we observe line formation from above or 
the side) and the angular dependence of the resonant scattering cross-section of the N harmonic: 

a^ « (1 + ^1) (1 - ^^lr-\ (4) 

If hb = 0, then the photon travels perpendicularly to B, and ajy is non-zero for all N, and a 
maximum for all iV > 1. If fis = !> then the photon travels along B, and a^ equals zero for all 
A^, except N = 1, for which case the cross-section is maximized. Radiative decay to the Landau 
ground state is an electric dipole transition, and the probability of photon emission in a particular 
direction is oc 1 + fi"^ (Daugherty & Ventura 1977); the emitted photon is thus most likely to be 
emitted in the direction of B. 

If we compare the line profiles in Figures 3-5, we see that they are broadened if the observer 
is oriented along B. This is because scattering will only occur if the energy of the photon, as 
calculated in the electron rest frame, matches the resonant energy, i.e. if -E^(l — [5^b) = -E-o where 
(5 = "^ , E^ \s the photon lab frame energy, and E^ is the line centroid energy for the first harmonic. 
Hence if ^b = 0, only those photons with exactly the resonant energy can scatter; if /i^ = 1, then 
there is a range of E^ such that scattering is possible, so that the scattering profile is broadened. 

If we compare Figures 3a and 4a, we see distinctive "shoulders" arising in the red and blue 
wings of the first harmonic in the 1-1 geometry. They appear most prominently when the observer 
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is oriented along the fi (p ^ 1), regardless of the angle ^ or the geometry (Isenberg et al. 1998b). 
Line shoulders are a predicted element of radiation transfer in strong fields (Wasserman & Salpeter), 
and many authors discuss their appearance (Alexander & Meszaros 1989, Nishimura & Ebisuzaki 
1992, Araya & Harding 1996, Isenberg et al. 1998b; of. Nishimura 1994, whose results contradict 
those of Alexander & Meszaros and Isenberg et al.). If oi, the ratio of the natural line width 
to the Doppler energy width for the first harmonic, -^ 0, then the energy of a scattered photon 
must lie within the regime of the optically thick line core. That photon will thus only escape after 
~ r^ scatters. If ai ^ 0, then with each scatter there is a probability ~ ai that a photon will 
be redistributed into the optically thin line wings, so that after ~ aj scatters, that photon will 
escape. Photons in the first harmonic core will, on average, scatter at least this many times for 
the 1-1 geometry, but not for the 1-0 geometry, for the values of A'^e considered in this paper. This 
is due to the fact that in the 1-0 geometry, photons which scatter in the line core may cross the 
source plane, whereupon they are lost from the calculation. The size of the shoulders will change 
depending upon the importance of Raman scattering in the line-forming region, since this dictates 
the number of photons spawned within the first harmonic line core. Their relative size in the red 
and blue wings is determined by the relationship between Tg ii and Tc; because we assume these 
temperatures to be equal, the equivalent widths of the two shoulders will be equal. 



3. Statistical Methodology 

In this work, we examine two general classes of models: the azimuthally symmetric model 
class, which has free parameters B, Nf., and ^u (^ = and is undefined); and a general model 
class which includes $ and as additional free parameters. If the GRB source is a neutron star 
with a simple dipole field, these two classes may be interpreted as representing line formation at 
the magnetic pole, or elsewhere. Hence we will refer to the azimuthally symmetric models as "polar 
cap" models throughout the remainder of this work. Because many physical processes are known 
to occur at the magnetic pole (e.g. accretion of gas from a companion star onto a neutron star), 
the association of simpler models with the magnetic pole has intuitive appeal. 

In this section, we summarize the methods of statistical inference which we use to fit models 
from each class to the data, to compare the best-fitting models from each class, and to place limits 
on values of the free parameters of the models which fit best overall. The reader will find fuller 
descriptions of these methods in Papers I and HI, and references therein. 
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3.1. Model Fitting 

We create Monte Carlo spectra at discrete points on a grid of values for the model parameters 
Bi2 and A'c,2i) and, if applicable, ^. We show these values in Table l}^ For each given set of 
parameter values, we determine the Compton temperature Tc (= T'eJj) W carrying out a series of 
Monte Carlo calculations for several different temperatures and determining at which temperature 
heating and cooling balance within the slab (see, e.g.. Freeman et al. 1992, Isenberg et al. 1998b). 
Hence, Tq is not a free parameter of the fit. 

Photons emerging from the line-forming region are binned. We generally use eight bins span- 
ning n = [0,1] (where /i = cos 6*; an exception made in fits to GRB870303 SI). As for 0, if ^ / 
and ^ / |, we may take advantage of symmetry across the {B, n)-plane to use eight bins spanning 
(j) = [0,7r]. If ^ = I, further symmetry across the plane perpendicular to the (5, n)-plane allows 
us to reduce the number of bins to four, spanning (p = [0,^]. 

We find that a spectrum is sufficiently accurate, i.e. the counts standard deviation in each 
output bin is sufficiently low, if ~ 10^ photons are injected into the line-forming region. We 
also find that the number of photons needed to accurately portray the low-energy line profiles is 
roughly an order of magnitude fewer than the number needed to accurately portray the high-energy 
continuum (E ^ 100 keV). We cannot ignore the high-energy continuum, because the Ginga GBD 
recorded the energy lost by an impinging photon (e.g. a 500 keV photon may have lost only 20 keV 
while passing through the GBD, leading to the recording of a count within the energy regime of 
interest). Thus we run the code with continuum photons sampled at energies ^ 100 keV, and attach 
to the resulting spectrum a separately created high-energy continuum (E ~ 100-1479 keV). Each 
spectrum is properly weighted so that overall spectrum is smoothly continuous at the matching 
point. 

For a given set of input parameters, many spectra will be generated, one for each binned value 
of (/U,(/>). In principle, we would assess the goodness-of-fit of each generated spectrum using the 
Poisson likelihood function C, and select that spectrum for which C is maximized. However, for 
reasons given below, we fall back upon the understanding that in the limit of a large number of 
counts n in a bin, we can use Pearson's x^ statistic, an approximation of L = log£: 

s' = ^ ("^--"-^ . (5) 

The fit extends over N data bins, and rrii and rii are the model amplitude and data in bin i, 
respectively. The best-fit parameters are those for which s^ is minimized. ^^ We set af = rui 



^^We can analytically shift the line profiles in E (or equivalently B12) by up to ± 10%, with no loss in accuracy, 
during fits. 

^■'We follow the notation of Larnpton, Margon, & Bowyer (1976), reserving the symbol x^ for a statistic which is 
explicitly sampled from the x^ distribution. 
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("model variance"), and hereafter denote the fitting statistic as s^. 



3.2. Model Comparison 

In Paper I, we describe both frequentist and Bayesian methods of model comparison, and 
apply both methods to determine the frequentist significance of, and the Bayesian odds favoring, 
a spectral model with exponentiated Gaussian absorption lines. The calculation of the Bayesian 
odds generally requires the analyst to numerically integrate the Poisson likelihood function L over 
parameter space. This can be computationally intensive when the number of parameters is large. 
However, we can determine the odds analytically if the shape of the likelihood surface in parameter 
space is similar to that of a multi-dimensional Gaussian (the Laplace approximation). In Paper 
I, we were able to reparametrize the exponentiated Gaussian line model so that we could use 
the Laplace approximation. Unfortunately, the function C{B, fi,N(>,^,(f>) does not have such the 
required Gaussian shape, nor can we reparametrize the model so that C will have this shape. Thus 
in this work we use only the frequentist method of model comparison. 

The frequentist comparison of two models, the null hypothesis Hq, and the alternative hy- 
pothesis Hi, is carried out by constructing a test statistic T, which is usually a function of the 
goodness-of-fit statistics for both models. There are two probability distribution functions (PDFs) 
which indicate the a priori probability that we would observe the value T, which are computed 
assuming the truth of Hq, and Hi, respectively. The test significance a is calculated by computing 
the tail integral of the Hq PDF from T to infinity. The resulting number represents the probability 
of selecting Hi when in fact Hq is correct; if a is sufficiently small, we reject Hq in favor of Hi. 
A common threshold for rejecting the null hypothesis is a < 0.05; in this work, we use the more 
conservative threshold a < 0.01. The selection of the threshold value is subjective; we use a more 
conservative value than is used normally because, e.g., we would demand that any increase in the 
quality of fit provided by general models be sufficiently great to justify considering the ramifications 
of line formation taking place outside the polar cap. 

The preferred model comparison statistic T would be the ratio of the maximum Poisson like- 
lihoods of Hi and Hq, but then simulations are required to estimate a. We can make an analytic 
estimate of a that is approximately correct using the x^ Maximum Likelihood Ratio (x^ MLR) 
test (Eadie et al. 1971, pp. 230-232). The test statistic is T = As^ = s^(i?o) - ^^(^i)- The PDF 
p{As'^\Hq) is then assumed to be the x^ distribution for AP = Pi — Pq degrees of freedom, where 
Pq and Pi are the number of free parameters in models Hq and Hi respectively. In order to use this 
test, the simpler model Hq must be nested within the more complicated model Hi, i.e. we arrive 
at Hq after setting the extra AP parameters of Hi to default values (often zero). 

For a single dataset, we directly compare fits of the three-parameter polar cap model, and 
five-parameter general model, to the data. When jointly fitting two (or more) datasets, the process 
of model selection becomes more complicated, as there is more than one model per class. For 
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instance, for fits to SI and S2, there are eight polar cap models, with the simplest specifying 
that (-Bsi) -^e,si) A*si) = (-632, -^6,32, A*S2)) and the next three simplest specifying that one of the 
parameters changes between the epochs of SI and S2, etc. (Analogously, there are 32 models in the 
general class.) We illustrate how we select models within a model class in Figure 6. We begin by 
comparing the simplest model (Mi) with all models that have a greater number of free parameters 
(M2-MN in Figure 6), computing the significance of the additional parameters of each alternative 
model (ai,2-ai,N)- We then adopt the simplest alternative model for which a^2i^Lj^ < 0.01 as our 
new null hypothesis Mi. (If two or more models with the same number of free parameters fulfills this 
criterion, we would select the one with the smallest s^, or equivalently, the smallest a;^2MLR-) We 
repeat the comparison process, until either no alternative model is adopted, or the most complex 
model is selected. 



3.3. Parameter Estimation 

In Paper I, we describe both frequentist and Bayesian methods of parameter estimation, and 
apply both methods to determine the confidence (frequentist) and credible (Bayesian) intervals on 
the free parameters of the best-fit exponentiated Gaussian absorption line models. If the surface 
defined by s^ in parameter space has paraboloidal shape, then the na confidence interval for 
an individual parameter is given by those values of that parameter for which s^ = Smmin + "^' 
where the values of all other free parameters are allowed to vary to new best-fit values. Otherwise, 
simulations are needed to determine the appropriate confidence intervals (Eadie et al., pp. 190-201). 
Since the s^ surfaces for fits in this work do not have the required paraboloidal shape, we limit 
ourselves to computing Bayesian credible intervals. 



We may determine a credible interval for a particular parameter x of model M, without 
reference to the other, "uninteresting," parameters, collectively denoted x', by marginalizing the 
Bayesian posterior function p{x,x'\D,I) over the space of parameters x': 

p{x\D,I)= dx'p{x,x'\D,I) ca dx'p{x,x'\I)p{D\x,x' ,1) . (6) 



D and / represent the data and background information about the experiment (such as the detector 
bandpass) respectively, while p{x, x'\I) is the prior probability (a quantitative statement of our state 
of knowledge about the relative probability of each possible value of x and x' before the data D are 
examined). p{D\x, x' , I) is simply the likelihood, C{x, x'), which is assumed to be oc exp(— s^/2). If 
the prior is assumed to be uniform, then p{x\D,I) is simply proportional to C{x,x'). The credible 
interval is then the range [xi,X2] such that 

!^'dxp{x\D,I) 



where z is the desired probability content (e.g. 0.683 for la bounds), and p{xi\D,I) = p{x2\D,I). 
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(Note that if the Hkehhood surface is not well-behaved and exhibits a number of modes, then the 
credible "interval" may actually consist of a number of intervals.) 

We note two items which affect the interpretation of these regions. First, our use of the Monte 
Carlo code is limited to iVc,2i ~ 10. The upper limit can potentially cut off highly probable regions of 
parameter space, affecting the computation of credible intervals. Second, it is not computationally 
feasible to create models using all \1' values when creating credible intervals for the parameters of 
the general model. Thus the quoted marginal distribution for each model parameter other than ^ 
is of the form p{x, ^ = ^best fitl^j I), and we do not compute the marginal distribution for ^ itself. 



4. Application of the Cyclotron Scattering Model to GRB870303 

4.1. The Data 

The Los Alamos/ISAS Gamma- Ray Burst Detector (GBD; Murakami et al. 1989) on Ginga 
detected GRB870303 at 16:23 UT on 3 March 1987. Figure 7 shows burst-mode time history data 
for the GBD Proportional Counter (PC), which covered 1.4-23.0 keV, and the GBD Scintillation 
Counter (SC), which covered 16.1-335 keV. The GBD continuously recorded burst-mode data at 
0.5-second intervals. These data were not stored in memory until a burst was detected, at which 
time the data from 16 seconds prior to the burst trigger until 48 seconds after the burst trigger were 
stored. The background-subtracted GRB870303 spectra (shown in Figure 8) exhibit lines in two 
time intervals, denoted Spectra 1 and 2 (SI and 82) by Graziani et al. (1992). SI is a 4 s spectrum 
with photon fluence 1.3 x 10~^ erg cm~^ in the bandpass 50-300 keV (estimated from the best-fit 
model of Paper I), that exhibits a saturated line at ~ 20 keV. S2 is a 9 s spectrum with fluence 4.5 
X 10~® erg cm~^ exhibiting harmonically spaced lines at ~ 20 and 40 keV. The midpoints of SI 
and S2 are separated by 22.5 s. 

The burst detector on Pioneer Venus Orbiter (PVO) also observed GRB870303, allowing the 
use photon time-of-arrival information to limit the possible source location to an annulus upon the 
sky. Yoshida (private communication, correcting Yoshida et al. 1989) reports this annulus to lie in 
the range 11.2° ^ ^inc ^ 57.6°. Since the shape and amplitude of a model counts spectrum that is 
derived from a given photon spectrum depends sensitively on ^inc, we treat this angle as a freely 
varying model parameter in Paper I. However, testing the cyclotron scattering hypothesis for more 
than one angle of incidence is not (currently) computationally feasible. Thus, we assume ^inc = 
37.7° in this work; this is the angle assumed by others who analyze GBD data (e.g. Murakami et 
al. 1988, Fenimore et al. 1988, Wang et al. 1989a, Graziani et al. 1992, Graziani et al. 1993). 
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4.2. Results: GRB870303 SI 



In Paper I, we describe fits to these data using exponentiated Gaussian absorption lines. We 
fit the data in PC bins 8-15 and SC bins 2-31 with a power-law model {-^ oc E~^'^^), onto which 
is multiplied one line with equivalent width We ~ 10 keV. For this model, s^ = 22.3 for 34 degrees 
of freedom (dof). The frequentist significance of this line is 1.2 x 10^^, and the Bayesian odds 
favoring the model with a line is 120-1. (These numbers differ from those given above, and in Paper 
I, in that Oiac is assumed to be 37.7°.) 

When we use the 1-0 geometry, the mode of the polar cap model places the observer along B, 
where the first harmonic equivalent width is maximized, while higher harmonic equivalent widths 
are minimized. For this model, s^ = 33.0 (33 dof), with B12 = 1.78, A^e,2i = 5, and 0.875 < // < 1. 
This fit compares unfavorably with the best phenomenological model fit (s^ = 22.3). We decrease 
A^ and examine ^ ~ 1, to maximize the first harmonic strength and to minimize the effect of the 
weak second harmonic. If we use a bin-width Afi = 2^, s^ = 26.7, with A'c,2i = 10 (the upper 
limit) (see Figure 9, and Tables 2 and 3). Smaller bin widths do not lead to further reduction in 
s^. We note that tests using values of A'e,2i above the upper limit indicate that s^ increases, due 
to increasingly apparent shoulders in the red and blue wings of the line. This fit indicates that we 
observe line formation from a privileged position above the polar cap if the line-formation region 
is illuminated isotropically. 

The difference As^ between this best-fit and that of Paper I is due to the fact that at the 
polar cap, we cannot create a first harmonic line with a sufficiently large equivalent width (VFe ~ 
10 keV). To determine the equivalent width of the best-fit line model, we model the line profile 
using three exponentiated Gaussians, two representing the shoulders and one representing the line 
core: 

^{E) = C{E) X n|=fexp(-Aexp[- ^^~J-''^' ]). (8) 

C{E) is the continuum photon flux at energy E, /?2 = /3core > 0, and /3i and [3^ < 0. The equivalent 
width of the line is We ~ 5 keV, a value that depends only weakly upon Nc in the vicinity of 
the mode. (This value is roughly 2a, or As^^ ~ 4, away from the best-fit phenomenological value 
We ~ 10 keV.) The equivalent width of the line core is ~ 7 keV, while each of the shoulders has 
equivalent width ~ — 1 keV. 

One can increase the equivalent width is to limit photon input into the line-forming region to a 
cone with axis parallel to n. This decreases the number of photons spawned into the first harmonic 
line core by greatly reducing the angle-averaged cross-section of the second and third harmonics. 
Such a model would be consistent with the hypothesis that photons are beamed during continuum 
formation (Share et al. 1986; Ho & Epstein 1989; Ho, Epstein, & Fenimore 1990; Dermer 1990). A 
typical beam cone opening angle for 20-40 keV continuum photons is ^cone ~ 45°. Such a model is 
also consistent with the hypothesis that the line-formation region is suspended far above the polar 
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cap (Dermer & Sturner 1991, Sturner & Dermer 1994).^^ In this case, only photons which travel 
along the magnetic polar axis will interact with the localized region of line formation (^cone ^0°). 
We compute spectra for these two test cases, assuming B12 = 1.76 and A'^e^2i = 0.6 (6% of the 
best-fit value for isotropic photon input). We find that both hypotheses are consistent with the 
data; in both cases, s^ = 23.0, and the line equivalent width is increased to « 10 keV. 

We find that the general cyclotron scattering model, with ^ and 4> as additional free param- 
eters, fits to the data better than the polar cap model, with s^^ = 23.6 (31 dof). The best-fit for 
this model occurs for ^ = -1; we may interpret this as meaning that line formation occurs at the 
neutron star magnetic equator. The observer is oriented along B (i.e., is looking along the surface 
of the slab); the line-of-sight opacity in the first harmonic is thus maximized (see Figures 3a and 
5b), and the equivalent width of the line can reach the desired value ~ 10 keV. 

Comparing the polar cap and general models, we find that the significance of the extra two 
parameters in the general model is o:~^2y[LR = 0.21; we select the polar cap model. It is however 
of great theoretical interest that models where lines form far from a polar cap can fit the data so 
well, a point to which we return in §5.1. 

For the 1-1 geometry, the polar cap and general models fit to the data with s^ = 37.0 (33 
dof) and 24.4 (31 dof), respectively. Polar cap model spectra do not fit well to the SI data because 
of the presence of either line shoulders (which have maximum size at /U = 1) or a strong second 
harmonic line (which has maximum strength at ^ = 0). Comparing the polar cap and general 
models, we find that the significance of the extra two parameters of the general model is ciy2^i^^ 
= 1.6 X 10~^; we select the general model (see Figure 10, and Tables 2 and 3). The SI data thus 
indicate strongly that line formation does not occur at the magnetic polar cap, within the context 
of a semi-infinite atmosphere. 



4.3. Results: GRB870303 S2 

In Paper I, we describe fits to these data using exponentiated Gaussian absorption lines. We 
fit the data in PC bins 7-15 and SC bins 2-31 with a power-law model that is cut off exponentially 
(ze '^ E^^'"^"^ X exp[ j^oyo ])• To model the lines, we multiply onto this continuum two lines with 
equivalent width We ~ 2.3 and 4.6 keV. For this model, s^^ = 33.6 (35 dof). The frequentist 
significance of these lines is 4 x 10^^, and the Bayesian odds favoring the model with lines is 
17-1. (Again, these numbers differ from those given above, and in Paper I, in that they are derived 
assuming Oi^c = 37.7°). 

For the 1-0 geometry, the polar cap and general models fit to the data with s^^ = 33.6 (33 
dof) and 31.7 (31 dof), respectively. We find that the relatively weak lines of S2 can form at any 



^■'Another, similar, hypothesis that cannot be directly tested with static models is that the line-forming region is 
driven off the star in a wind by radiation pressure (Miller et al. 1991, 1992; Isenberg et al. 1998a). 
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location on a magnetized neutron star, i.e. for each ^, there exists some values of /i and (p such 
that the data are fit well with the cyclotron scattering model. The significance of the two extra 
parameters of the general model is a^2MLR = 0.39; we select the polar cap model (see Figure 11, 
and Tables 2 and 4). 

For the 1-1 geometry, the polar cap and general models fit to the data with s^ = 33.9 (33 dof) 
and 32.9 (31 dof), respectively. Again, as is the case for the 1-0 geometry, we find that the S2 lines 
can form anywhere on a magnetized neutron star. The significance of the two extra parameters of 
the general model is a^2^i^^^ = 0.61; we select the polar cap model (see Figure 12, and Tables 2 
and 4). We note that ^ ^ for this geometry, unlike for the 1-0 geometry; this is directly related 
to the presence of shoulders in spectra as /i ^ 1. 



4.4. Results: the Combined (81-^82) Data 

In Paper I, we describe fits to combined (81+82) data using exponentiated Gaussian absorption 
lines. We select a line model parameterized by the first harmonic energy (Ec ~ 21.6 keV), two 
first harmonic equivalent widths (VFe,i,si ~ 10.3 keV and VFe,i,S2 ~ 2.5 keV), and one second 
harmonic equivalent width (VFe,2,si = 1^E,2,S2 ~ 3.1 keV). For this model, s^^ = 58.0 (68 dof). The 
frequentist significance of the lines, evaluated jointly, is 3.1 x 10^^, and the Bayesian odds favoring 
the model with lines is 8080-1. (Again, these numbers differ from those given above, and in Paper 
I, in that they are derived assuming 6i^c = 37.7°). 

As noted in §3.2, there are eight polar cap models for joint fits to two datasets, which have 
a minimum of three, and maximum of six, free parameters. After selecting those models with 
the lowest values of s^ for each possible number of free parameters, we find s^ = 70.8 (69 dof), 
66.5 (68 dof), 64.4 (67 dof), and 60.3 (66 dof), respectively for the 1-0 geometry (Table 5). In no 
case is As^ sufficiently large so that we may reject the null hypothesis; the significance of the 3 
additional parameters of the most complex model is 0^2^11^^ = 0.015. The additional reduction in 
s^ that would be necessary for a^2MLR to reach our selection criterion of 0.01 is ~ 0.9. Because the 
variation of s^^ between Monte Carlo manifestations of input models is ^ 1, there is only a small 
possibility that the 6-parameter model would be selected more often that the 3-parameter model 
if we were create and fit an infinite number of Monte Carlo spectra. Thus we select the simplest 
3-parameter polar cap model (see Tables 2 and 6, and Figure 13; the credible regions, which we do 
not present, are similar to those in Figure 11). 

Fits of the 32 general models (with minimum five, and maximum ten, free parameters) to the 
combined (81+82) data for the 1-0 geometry yield s^^ = 69.5 (67 dof) for the simplest 5-parameter 
model, and 59.6 (66 dof) for the 6-parameter model for which /xsi 7^ AiS2 (Table 5). The significance 
of the additional parameter is OyZMLR = ^-"^ ^ 10"^; we select the 6-parameter model. No more 
complicated model is favored over this model (see Tables 2 and 6, and Figure 14; in this figure, we 
show only the credible regions for the parameter subset [/Usi,;Us2])- 
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Comparing best-fit polar cap and general models, we find that As^ = 11.2 for AP = 3. The 
significance of the additional parameters is CKy^MLR = 0.01. We conclude that on the basis of the 
current evidence, we cannot choose between these models. 

The application of the polar cap models to the combined (S1+S2) data for the 1-1 geometry 
yields best-fits of s^ = 73.7 (69 dof) for the simplest model, through s^ = 71.3 (66 dof) for the 
most complex model; a^2-^i^Yi = 0-49 (Table 5). We select the simplest model 3-parameter model 
(see Tables 2 and 7, and Figure 15; the credible regions, which we do not present, are similar to 
those in Figure 11). 

The application of the general class of models yields s^ = 71.6 (67 dof), 65.0 (66 dof; fisi / 
fiS2), and 59.9 (65 dof; /isi / fJ'S2 and (psi / 032), respectively; fitting models with eight or more free 
parameters offers little further reduction in s^. The significance of the first additional parameter is 
Q^x^MLR = 0.01; we cannot choose between the 5- and 6-parameter models. We describe the result 
of choosing each of these models in turn: 

1. If we choose the 5-parameter model, we would compare it to the 7-parameter model and find 
Q^y^MLR = 2.9 X 10"'^. We would select the 7-parameter model. We would compare this 
model to the best-fit 3-parameter polar cap model and determine that As^^ = 13.8 for AP = 
4; a^2j^Lj^ = 8.0 x 10^^. The amount by which As'^^ surpasses the a-^ajviLR = 0.01 criterion 
is ~ 0.5. This is insufficient to conclude that the data select the more complex model, and 
we would conclude that we cannot select between models. 

2. If we choose the 6-parameter model, we would compare it to the 7-parameter model and find 
a^2MLR = 0.024. We would select the 6-parameter model. We would compare this model to 
the best-fit 3-parameter polar cap model and determine that As^ = 8.7 for AP = 3; a^2y[i^Yi 
= 0.034. We would select the 3-parameter polar cap model. 

Thus, again, we find that we cannot select between the simplest 3-parameter polar cap model and a 
7-parameter general model on the basis of the current evidence. (Best-fit values, credible intervals, 
etc., for the 7-parameter general model may be found in Tables 2 and 7, and Figure 16; note that 
in this figure, we show only the credible regions for the parameter subset [/isi) A*S2 5 '/'sij ';^S2]-) 



4.5. Limits on Neutron Star Rotation Period 

The combined (S1+S2) data indicate that if ^ = f, the best-fit values of fi and/or 4> change 
during the 22.5 s between SI and S2 for both the 1-0 and 1-1 geometries. In other words, the 
orientation of the observer relative to the line-forming region changes with time. The simplest way 
to interpret this result is to invoke neutron star rotation (see, e.g.. Lamb, Wang, & Wasserman 
1992). 

For the 1-0 geometry, the best-fit model has fi changing as a function of time, but not (f> 
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(Table 6). Because of the symmetries which exist for the particular case of ^ = ^, </) G [~f ) f ] 
or G [x'T"l" However, if ^ or vr, then (\) must change as ii changes, violating the original 
model assumption that (/> does not change. This condition leads us to posit a model of the rotating 
neutron star in which: (1) the rotation axis is perpendicular to the observer's line-of-sight; (2) the 
star may rotate in either the clockwise or counter-clockwise directions; (3) the magnetic axis lies 
within the neutron star-observer plane; and (4) line formation is localized and occurs only where 
the magnetic equator intersects the neutron star-observer plane. This model is shown pictorially 
in Figure 17. We assume that the neutron star rotates less than once between SI and S2, since 
otherwise we could expect to see evidence for line(s) during the time period separating the two 
spectra. 

The rotation period is given by 

*-t = \a '^'^ a I X 22.5s, (9) 

where 6*31 = cos~^(/Usi) and 6*32 = cos~^(^s2)- Generally, there are two possible rotation periods 
that we may derive, depending upon the direction of the star's rotation. However, for ^ = ^, there 
are four possible times (which we hereafter denote as trot,n), because of the symmetry between 

= and = vr (see Figure 17). 

To define credible intervals for the rotation periods, we repeatedly sample values of the cosines 
/isi and /is2 from the posterior distribution p(/xsi,MS2|^ = I, = or vr, D) and determine for 
each sampled cosine pair the possible rotation periods trot,n- (This distribution is slightly different 
from that shown in Figure 14, which includes marginalization over 0.) We determine probability 
distributions from 10^ values of trot.n, and use these distributions to estimate the 1, 2, and Zo 
credible intervals given in Table 8. 

For the 1-1 geometry, both [i and (/> change as a function of time, which greatly complicates the 
derivation of rotation periods. In particular, the orientations of the magnetic and rotation axes, H 
and i?, can be arbitrary. In the Appendix, we describe a Bayesian method which we use to derive 
rotation periods once we sample values of /i and ^ for SI and S2. We determine the probability 
distributions from 10^ values of trot, id and use these distributions to estimate the 1, 2, and 3cj 
credible intervals given in Table 8. 



5. Discussion 

In this paper, we demonstrate that Monte Carlo models of cyclotron scattering in the strong 
magnetic field {B ~ 10^^ G) of a galactic neutron star can successfully account for the positions and 
strengths of the fines exfiibited at f» 20 keV in GRB870303 SI and « 20 and 40 keV in GRB870303 
S2. Our results are robust to cfianges in slab geometry and magnetic field orientation. Given that 
physically rigorous models of line formation within the cosmological burst environment do not yet 
exist (Stanek et al. 1993 and Ulmer & Goodman 1995, e.g., invoke femtolensing) , our results, when 
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paired with the successful fits of cyclotron scattering models to the data of GRB880205 by Wang 
et al. (1989a) and Freeman et al. (1992), provide strong support for the hypothesis that that some 
(though not all) GRBs are galactic in origin. 



5.1. Galactic Source Population 

5.1.1. Static Line Formation at the Polar Cap 

We model static line formation at the magnetic polar cap of a neutron star with a simple 
dipole field by setting B parallel to the slab normal n, i.e. by setting ^ = 0. We find that while we 
can easily fit this model to the data of GRB870303 S2, acceptable fits to the data of GRB870303 
SI can be made only with difficulty (1-0 geometry), or cannot be made at all (1-1 geometry). For 
the 1-0 geometry, we conclude that either the observer is oriented directly along B (p -^ 1), or 
that the line-formation region is levitating far above the polar cap, where it presents a small solid 
angle to continuum photons, with the result that the equivalent widths of the second and third 
harmonics are reduced relative to the width of the first harmonic. The latter conclusion is consistent 
with the suggestion by Dermer &; Sturner (1991), and Sturner & Dermer (1994), that a scattering 
atmosphere with a geometry similar to the 1-0 may be present in both accretion-powered pulsars 
and GRBs, with layers that are optically thick to line scattering but thin to continuum scattering 
forming within a few stellar radii of the neutron star center. For the 1-1 geometry, we find that 
prominent emission-like shoulders on either side of the first harmonic prevent a good fit to the SI 
data if the observer is oriented along B (/U = 1); by decreasing /i, we can reduce the magnitude of 
the shoulders, but second and third harmonics form, preventing an acceptable fit to these data.^^ 

The results of fits to the SI data argue against static line formation at the magnetic polar cap. 
If line formation does indeed occur above the polar cap in an outflowing plasma, the neutron star 
probably resides within the galactic halo at a distance of ^ 50 kpc. Two lines of reasoning point 
to this conclusion. First, as pointed out by Lamb et al. (1990), a static polar cap line formation 
region will be disrupted on timescales ~ 10^® s unless the the burst flux is below the magnetic 
Eddington limit. At the non-magnetic Eddington limit, the radiative pressure exerted by both line 
and continuum photons upon electrons in the line-forming layer balances the gravitational force 
upon the protons which are electrostatically coupled to the electrons; the magnetic Eddington limit 
is determined by replacing the Thomson cross-section of the electron with the resonant cyclotron 
scattering cross-section of the first harmonic. This lowers the luminosity limit from ~ 10^^ erg s^^ 
to ~ 10^^ erg s~^. The upper limit on the distance to the burst source with an electron-proton 
line-forming layer is thus ~ 100 pc. 



^^However, we note that shoulders are not always inconsistent with observed GRB line profiles; Freeman et 
al. (1992) show that the 1-1 polar cap model, with shoulders, fits to the data of GRB880205 better than a 1-0 polar 
cap model without shoulders; s^ falls from 43.9 to 36.6. 
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Second, Loredo &; Wasserman (1998b) determine that the data of the 3B catalog (Meegan 
et al. 1996) is consistent with the hypothesis that there is a component of the source population 
residing within the galaxy, with that component either being comprised of dim local halo sources 
at distances ^ 1 kpc, or luminous halo sources at distances ^ 50 kpc. They make this conclusion 
by using Bayesian methods (developed in Loredo & Wasserman 1995) to compare a model in 
which sources residing in a Bahcall-Soneira halo with core size 2 kpc (Bahcall & Soneira 1980) 
are mixed with cosmological GRB sources, with the best-fit cosmological GRB source model from 
Loredo & Wasserman (1998a). Each class is assumed to contain standard candle bursters, with the 
further assumptions for the cosmological bursters that the comoving burst rate is homogeneous and 
independent of redshift. (Loredo & Wasserman 1998a demonstrate that this simple cosmological 
model fits the data as well as more complicated models with either an inhomogeneous burst rate 
or a power-law burst luminosity function.) The best fit of the halo model is better than that of the 
purely cosmological model, for both dim local and luminous halo bursters. The best-fit fractions of 
dim local halo sources are 0.59 (64 ms 3B catalog data) and 0.36 (1024 ms data), and the bursters 
are limited to distances ^ 1 kpc; the respective fractions for luminous halo sources are 0.11 and 
0.073, with inferred distances ^ 50 kpc. However, the 3B data lack the power to tightly constrain 
the model parameters, and thus the credible intervals for each fraction include zero at a significance 
level ^ 2o", or 95.5%. (Despite the lack of constraining power, the data do not support the bursters 
residing at distances intermediate between 1 and 50 kpc.) Because of this lack of constraining 
power, the Bayesian odds favoring the galactic component models are not large: 2.5 and 6.7 (64 
and 1024 ms) for dim local halo sources, and 0.45 and 0.25 for luminous halo sources. An odds of 
> 10-20 would be considered strong evidence in favor of an alternative model (see the review by 
Kass & Raftery 1995 and references therein). 

We note that if lines exhibited by accretion-powered pulsar (APP) spectra are also formed 
via cyclotron scattering, the lack of shoulders in the 12 known APPs with lines (Mihara 1995) 
places a lower limit on column densities that is much larger than the column densities determined 
in our analyses, within the context of the 1-1 geometry. This result is shown by Isenberg et 
al. (1998b), who consider cylindrical line-formation regions with photon injection along the axis, 
which represent the canonical model of the emission region of accretion-powered pulsars. They 
show that if Ne^2i ^ 1000, the shoulders disappear. This implies that even if the 1-1 geometry is 
an acceptable approximation for both APPs and GRBs, the mechanisms underlying line formation 
must be significantly different. 



5.1.2. Static Line Formation at the Magnetic Equator 

The most theoretically intriguing result of our analyses is the possibility that line-formation 
may occur away from the magnetic polar cap (i.e. for ^ 7^ 0). The SI data provide strong evidence 
that, for the 1-1 geometry, line formation does not occur at the polar cap, with the best-fit location 
being the magnetic equator (^ = ^). The combined (S1-I-S2) data provide marginal evidence, for 
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both geometries, favoring equatorial line formation. Line formation at a magnetic equator may be 
static even if the burst luminosity greatly exceeds the critical Eddington luminosity in a strong 
magnetic field, since the closed magnetic field lines there may effectively trap the plasma (see, e.g., 
Zheleznyakov & Serber 1994, 1995). Thus the assumption of line formation in static layers may be 
acceptable even if the progenitor neutron star resides in the Galactic halo at distances ^ 50 kpc. 



5.1.3. Line Formation in an Outflow 

If line formation occurs at the magnetic polar cap of a luminous halo burster, the line formation 
region will initially flow outwards along B. In an outflow, the variation of the magnetic field and 
plasma velocity with altitude will tend to broaden the lines. Thus the relevant question is whether 
we can still observe narrow lines in the spectra of distant halo GRBs. The answer appears to 
be yes. Miller et al. (1991,1992) calculated the properties of the second and third harmonics, 
approximating Raman scattering with cyclotron absorption (see §5.2). They showed that narrow 
lines can be formed at these harmonics, and they successfully fit outflow absorption models to the 
second harmonic of GRB880205. Chernenko &: Mitrofanov (1995) calculate the properties of the 
first harmonic line in an outflow, but they assume the line to be formed via absorption; they find that 
narrow first harmonic lines are possible. Isenberg et al. (1998a) show that the cyclotron scattering 
model can create narrow lines at all harmonics. They treat this problem using a variant of the 
Monte Carlo code applied in this work. They assume a hot spot at the magnetic pole of a neutron 
star, and apply the radiation force calculation of Mitrofanov & Tsygan (1982) to determine the 
plasma outflow velocity. They conclude that cyclotron scattering lines can form within a relativistic 
outflow with properties similar to those of GRB870303 and GRB880205, provided the hot spot is 
a small fraction of the stellar surface (^hot ^ O.li^Ns)- 

A natural consequence of the fact that the line-forming layer may move along field lines is 
that the height of the line-forming layer may change with time, as, for example, the burst intensity 
or hot spot radius changes. A change in height would be accompanied by a change in inferred 
magnetic field strength (i.e. by an inferred shift in line-centroid energy). Yoshida et al. (1992) 
report that the magnetic field strength shows a declining trend within the 9 s interval in which the 
line candidates of GRB880205 have maximum significance, which indicates that the line-forming 
region is flowing outwards from the neutron star. We find that the combined (S1+S2) data are 
nearly consistent with the hypothesis that the line-forming region moves inward during the ~ 20 s 
between line epochs. Within the context of the most complex 6-parameter polar cap model, which 
the data almost favor over the simplest 3-parameter model, B12 changes from 1.76 to 1.96 from SI 
to S2, indicating that S2 line formation might occur closer to the polar cap than SI line formation. 
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5.2. Cyclotron Scattering versus Cyclotron Absorption 

Many authors (e.g. Fenimore et al. 1988, Graziani et al. 1992) fit data with the cyclotron 
absorption model, motivated by its computational simplicity and by the fact that approximating 
cyclotron Raman scattering by cyclotron absorption is approximately valid for the second and 
higher harmonics. The cyclotron absorption model spectrum is 



CA{E) = C{E)exp 



m 
Y^AnGniE) 
n=l 



(10) 



where C is the continuum spectrum, m is the number of harmonics, and 



«»<^' = i-K-^^l <"' 



4-1, 2kT II u'- 

is the absorption line profile for the n harmonic. The line widths AEn are given hy E„- ' — — — 



where /x' is the cosine of the angle between the line of sight and B. The line amplitudes An are 
given by A*"]" ^On, where A^gn i^ the column density along the line of sight, and 

5.5 X 10-20 /n2 BY-\^ , , ,2^,, . /2^n-l ,„.. _-2 



is the absorption coefficient of the n**^ harmonic. Because modeling the first harmonic line with 
cyclotron absorption is not valid, we must allow different temperatures and column densities when 
fitting to the first, and higher, harmonics. If, for example, we fit two harmonic lines to the data, 
the relationships between the cyclotron absorption model parameters and the physical parameters 
of the line-forming region are Ei oc B, Ai oc Nl^Kl+fi'"^) and A2 oc i?Af^°|(l— ^'^), and AEi = 



£^n/ '^■ii-T and A^2 = 2^n/ "-"-T • iV'T and '=■"■7 , however, have no direct physical 
meaning. 

We may use the cyclotron absorption model to describe the line-formation region if at least 
both the second and third harmonics are strong, i.e. the data request that each of these harmonics 
be fit with a line shape parameterized by centroid energy, equivalent width, and full-width (see 
Paper I for a description of this particular parameterization) . Only then is the number of required 
fit parameters, five (we assume £^3 = §-£2)) larger than the number of physical parameters needed to 
describe the region {B, Tg y 2 = ^e,|l,3) -^1°! = ^1°3> /"')• Otherwise, such as when only the first and 
second harmonics are observed, radiative transfer calculations are required. No simple description 
can be used to explain the first harmonic line, the appearance of which depends critically on the 
outcome of the multiple resonant scatters required in order for individual photons to escape, as 
well as on the introduction of spawned photons at energies near that of the first harmonic. 

The cyclotron scattering model can describe the line-forming region using fewer free parameters 
than the cyclotron absorption model (because of the relationship between Tg m and B in the scat- 
tering model), and can use the first harmonic data to help forge the description. However, we have 
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found that moderate resolution GRB data may require an even smaller number of free parameters 
than the minimum three of the cyclotron scattering model. In Table 9, we list the number of model 
parameters required in fits to the GRB870303 data with the exponentiated Gaussian absorption 
line model (Paper I), and the cyclotron absorption and scattering models. Each harmonic in these 
data can be adequately modeled with a saturated line shape parameterized by centroid energy and 
equivalent width; the full-width is set to be proportional to the equivalent width. The number of 
free parameters is then A^harm + 1) where A'^fiarm is the number of observed harmonics (A^harm line 
width parameters and one harmonic line energy). This is an upper limit. For instance, the S2 data 
is adequately fit using two, rather than three, parameters (by setting We,2=2We,i), because the 
second harmonic data lack constraining power. 

Lamb (1992) argues that physically-based models which have parameters not required by the 
data will (1) adequately fit the data, but (2) have reduced diagnostic power because of a lack of 
constraint on the individual model parameters. This is certainly true for the cyclotron absorption 
model. For both GRB870303 (Graziani et al. 1992) and GRB880205 (Fenimore et al. 1988), the 
cyclotron absorption model determines the magnetic field well (since it is a function of only Ei), 
but n' is undetermined and kTf,u and Nl°^ are therefore poorly constrained. However, we determine 
in this work that photon spawning and the presence of line shoulders give the cyclotron scattering 
model greater diagnostic power than we would have predicted. (Note that such the absorption 
model fits adequately because large shoulders are not observed in GRB870303 and GRB880205; 
if large shoulders were present, no model using absorption-like profiles would adequately fit the 
data.) Spawned photons fill the first harmonic in the polar cap model, decreasing its ability to 
adequately fit the SI data in the 1-0 geometry; the addition of shoulders in the 1-1 geometry makes 
no adequate fit possible. While polar cap models adequately fit the S2 data, the line shoulders 
which appear in the 1-1 geometry greatly reduce the range of fi consistent with the data (cf. the 
1-0 result, for which there is no constraint on n at the 3a limit). 

Fits to the combined (S1-I-S2) data in Paper I clearly indicate that the data prefer a four- 
parameter model in which We,! changes between SI and S2. A cyclotron scattering model having 
different /U and/or A'c values for SI and S2, and at least four parameters overall, would thus be 
expected to provide the best fit to the data. Fits of the general model to the joint data fulfill 
that expectation, with six and seven parameter models chosen. Surprisingly, however, the simplest 
three-parameter polar cap models adequately fit the data for both the 1-0 and 1-1 geometries; 
spawned photons and line shoulders prevent the data from selecting more complex models over the 
simplest model, even if up to six parameters are allowed to fioat freely (Table 5). 



5.3. Neutron Star Rotation 

The use of cyclotron scattering models allows us to place rigorous limits upon the rotation 
period of the neutron star source of GRB870303. Previously, Graziani et al. (1992) interpreted the 
results of fits to GRB870303 SI and S2 to qualitatively place limits on the rotation period. They 
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use the observation of one harmonic during SI and two harmonics during S2 to suggest that hue 
formation is observed along B for SI and perpendicular to B for S2. Invoking both localized line 
formation and line formation along an equatorial arc, they derive a rotation period 45 s ^ irot ^ 
180 s. Lamb et al. (1992) invokes similar arguments to explain the changes in line-centroid energy, 
strength, and width of the line candidate exhibited by the HEAO-1 A4 data of GRB780325 (Hueter 
1987), and to predict 40 s ^ trot ^ 80 s. 

Examining the real-time data of GRB870303, Yoshida et al. (1989) find an additional peak in 
the PC time history beyond the two seen in burst mode (Figure 7). The peaks have periodicity ~ 
30 s. They conclude that the data in the third peak can be adequately fit an ~ 10 keV thermal 
bremmstrahlung spectrum, an ~ 1.6 keV blackbody spectrum, or an ~ 0.3 keV thermal cyclotron 
spectrum are consistent with these data. No spectral features are apparent, but the spectrum 
is perhaps too soft for lines to be unambiguously detected. Such a periodicity is perhaps barely 
consistent with trot, 3 in the 1-0 geometry, but is consistent with both trot, 3 and trot, 4 in the 1-1 
geometry. If we hypothesize that the neutron star completes more than half a rotation between 
SI and 82, it follows that a localized region of line formation disappears behind the stellar disk 
between the two epochs. Graziani et al. (1993) consider a model in which lines appear at time ti, 
disappear at time t2, and are constant in between. They find that the uncertainties at^ and 0"t2 in 
the times at which the lines appear and disappear are large, and they cannot exclude the possibility 
that lines are present throughout GRB870303. 

Rotation periods ^ 22.5 s are longer than those which would be expected for a neutron star 
born with rotation period 1 s and which has undergone magnetic braking for a Hubble time (.Pmax ~ 
= 5 s). Mass-transfer from companions slows some accretion-powered X-ray pulsars, which have 
rotation periods up to ~ 800 s (see, e.g., Nagase 1989). Similarly, accretion of matter onto high- 
velocity neutron stars on the way to or within a Galactic corona may act to both slow a neutron 
star and provide fuel for GRBs. Repeated GRBs (Graziani et al. 1998) may in turn cause losses in 
angular momentum, slowing the neutron star, if mass is ejected from the source. 

The authors would like to thank the referee, David Band, for his careful reading of the text 
and his many helpful comments. We also wish to thank Mike Isenberg, Carlo Graziani, and Tom 
Loredo for helpful discussions. This work was supported in part by NASA Graduate Traineeship 
NGT-50778 and NASA Grant NAGW-830. 



A. General Derivation of the Neutron Star Rotation Period 

The orientation of the observer relative to the slab where line formation takes place is described 
by a polar cosine, fi, and azimuth angle, 4> (see Figure 1). If both fi and (p change with time, we 
cannot use the simple model described in §4.5, in which cp is constant in time, to place limits on 
the rotation period(s) of the underlying neutron star. Here, we outline a general Bayesian method 
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for deriving credible intervals for the rotation period(s). 

We first write out a marginalized posterior probability distribution for fi and (p for the best- 
fit model, adding as free parameters the magnetic field axis orientation at the epoch of the first 
spectrum, Hi, and the rotation axis, R: 



I I dBdN,p{B\I)p{N,\I) X (Al) 

>C(//i, M2, 01, </'2, i^i, -R, S, A^c, ^ = f ) 



p{D\I) 

(Here, we assume that ^ = ^; we do not marginalize over this parameter.) We determine the 
rotation periods trot.n by sampling from this posterior distribution. (There are four possible rotation 
periods that can be derived, for reasons given in §4.5 and illustrated in Figure 17.) However, we 
must deal with the Bayesian prior, the first factor on the right-hand side of eq. (A2). We expand 
it: 

p{lJ'i,l^2,(/)i,(l)2,Hi,R\I) = p{R\I)p{Hi\R,I) ■ ■■p{(/)2\(l)i,fi2,l^i,Hi,R,I) 

= p{R\I)p{Hi\I)p{fii\I)p{p2\f^i,Hi,R,I)x (A2) 

p{(j)i\fii,Hi, I)p{4>2\p.2, Hi,R, I). 

We randomly sample pi, p2, Hi, and R from uniform distributions. Because of polar symmetry, the 
values fii and p2 define circles centered on the visible disk of the neutron star where line formation 
may take place. Given Hi, we can define the magnetic equator, which either intersects the line- 
formation circle twice, or not at all. Since we assume ^ = ^, the vector pointing out of the star at 
the intersection points is simply hi (Figure 1). At the intersection points, O ■ hi = /Ui, where O is 
the vector pointing to the observer. If we define O as (0,-1,0), pi = —ni^y. The only unknown is 
then ni^x, which we can solve for numerically. Given the location of line-formation, we project O 
onto the local surface tangent (i.e. onto a vector which is perpendicular to hi). The dot product 
of this projected vector with the vector —Hi (the direction of the field at the magnetic equator) is 
the cosine of (pi . 

We may determine h2 by using the information that R- h is constant and O ■ h2 = H2- Either 
two solutions to these equations exist, or none at all. Given a solution, we next determine the 
possible local field direction — i^2- The vectors R, hi, and n2 allow us to use the law of cosines 
to determine P, the angle through which the neutron star rotates between the two epochs. We 
transform the Cartesian coordinate system such that R' = (0,0,1), rotate the system through angle 
/3, and undo the transformation. We then solve for (/>2, as above. 

The likelihood C does not depend upon Hi or R, so the integrated expression in eq. A2 is 
proportional to the four-dimensional posterior probability distribution p{pi, fi2, 01, (f'2\D, I,^ = ^). 
For the particular case of interest in this paper, this distribution is known: two dimensional "slices" 
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of it are shown in Figure 16. We use this distribution as a rejection function: we assess the relative 
probabiUty of our solution with respect to the best-fit parameters of the model, sample a random 



Psolutio 

rotation period t^ot.n- 



number r, and accept the solution if r < Hsoiution ^ jf ^^ accept the solution, we then solve for the 

Pbcst fit 
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Fig. 1. — An illustration showing the coordinate system used in this work. The infinite plane- 
parallel slab is threaded by a magnetic field B oriented at an angle ^ relative to the slab normal 
h. The angles 9 (or angle cosine fi = cos 9) and (j) specify the direction of escaping photons relative 
to B and h. 

Fig. 2. — An illustration showing the two plane-parallel slab geometries which we examine in 
this work, denoted "1-0" and "1-1." The numbers represent the relative electron column densities 
above and below the continuum photon source plane. Injected photons travel through the slab, and 
emerge from one of the plane-parallel faces of the slab. The observer is situated above the source 
plane. 

Fig. 3. — Sample cyclotron scattering line spectra generated using the GRB870303 S2 continuum 
spectrum, with B12 = 1.7 and Ne^2i = 0.6, for the 1-0 geometry and ^ = 0. In spectrum (a), the 
observer is oriented directly above the slab, while in spectrum (b), the observer looks along the slab 
surface. The line-of-sight column density and the scattering cross-section a^ oc (1-|-/xb)(1— A*b)^~^' 
where ^b is the cosine of the angle between observer and field, dictate the strength of the N^^ 
harmonic line. 

Fig. 4. — Sample cyclotron scattering line spectra generated under the same conditions as Figure 
3, but for the 1-1 geometry. Shoulders appear at the first harmonic in spectrum (a) for reasons 
described in §2.3. 

Fig. 5. — Sample cyclotron scattering line spectra generated for similar conditions as Figure 3, but 
for ^ = ^. As indicated by the diagram at left, in spectrum (a), the observer is oriented directly 
above the slab; in spectrum (b), the observer looks along the slab surface and along the magnetic 
field; and in spectrum (c), the observer looks along the slab surface, perpendicular to the magnetic 
field. 

Fig. 6. — A flow chart showing the method which we use to select the best fitting polar cap (^ = 0) 
and general (^ 7^ 0) cyclotron scattering models, when jointly fitting two (or more) datasets. We 
then directly compare the selected polar cap and general models with the x^ Maximum Likelihood 
Ratio model comparison test. See §3.2 for details. 

Fig. 7. — Ginga Proportional Counter (PC; top) and Scintillation Counter (SC; bottom) time 
histories of GRB870303. The PC data is presented in 1 s bins, while the SC data is presented in 0.5 
s bins. The burst triggered the recording of Ginga GBD burst-mode data at ~ 16 s; the preceding 
16 s of burst-mode data, in memory at the time of trigger, were recorded and not overwritten. 
Overall, 64 s of burst-mode data were recorded. Epochs SI (4 s) and S2 (9 s) are shown; the 
midpoints of SI and S2 are separated by 22.5 s. 
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Fig. 8. — Ginga GBD count-rate spectra for intervals SI and S2 of GRB870303, normalized by 
energy-loss bin width. SI exhibits a single line at « 20 keV, while the spectrum S2 exhibits 
harmonically spaced lines at ~ 20 keV and 40 keV. 

Fig. 9. — Left: The data and predicted counts (upper panel), photon spectrum (middle panel), and 
residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to 
the data of GRB870303 SI for the 1-0 geometry. Right: Two-dimensional la, 2a, and 3a Bayesian 
credible regions for this fit, as a function of the stated parameters. 

Fig. 10. — Left: The data and predicted counts (upper panel), photon spectrum (middle panel), 
and residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit 
to the data of GRB870303 SI for the 1-1 geometry. Center and Right: Two-dimensional la, 2a, 
and 3(7 Bayesian credible regions for the fit. 

Fig. 11. — Left: The data and predicted counts (upper panel), photon spectrum (middle panel), 
and residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to 
the data of GRB870303 S2 for the 1-0 geometry. Right: Two-dimensional la, 2a, and 3a Bayesian 
credible regions for this fit, as a function of the stated parameters. 

Fig. 12. — Left: The data and predicted counts (upper panel), photon spectrum (middle panel), 
and residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to 
the data of GRB870303 S2 for the 1-1 geometry. Right: Two-dimensional la, 2a, and 3a Bayesian 
credible regions for this fit, as a function of the stated parameters. 

Fig. 13. — The data and predicted counts (upper panel), photon spectrum (middle panel), and 
residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to the 
data of GRB870303 Sl-hS2 for the 1-0 geometry and ^ = 0. Left: The fit to SI. Right: The fit to 
S2. 

Fig. 14. — The data and predicted counts (upper panel), photon spectrum (middle panel), and 
residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to the 
data of GRB870303 S1-FS2 for the 1-0 geometry and ^ = f . Left: The fit to SI. Center: The fit 
to S2. Right: Two-dimensional la, 2a, and 3a Bayesian credible regions for the parameter subset 
(MSl,AiS2)- 

Fig. 15. — The data and predicted counts (upper panel), photon spectrum (middle panel), and 
residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to the 
data of GRB870303 S1+S2 for the 1-1 geometry and ^ = 0. Left: The fit to SI. Right: The fit to 
S2. 
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Fig. 16. — The data and predicted counts (upper panel), photon spectrum (middle panel), and 
residuals of the fit in units of a (lower panel), for the best cyclotron scattering line model fit to the 
data of GRB870303 S1+S2 for the 1-1 geometry and ^ = f . Left: The fit to SI. Left Center: The 
fit to S2. Right Center and Right: Two-dimensional Icr, 2(T, and 3cr Bayesian credible regions for 
the parameter subset (/Usi,;Us2, i;^si) 0S2)- 



Fig. 17. — A model of neutron star rotation for the 1-0 geometry, for '^ = '^. The rotation axis (7 is 
perpendicular to the observer's line-of-sight, and the star rotates in either the clockwise or counter- 
clockwise directions. The magnetic field axis lies within the neutron star-observer plane. Line 
formation is localized and occurs where the magnetic equator intercepts the neutron star-observer 
plane. Because of the symmetry between = and tt, there are four angles through which the 
neutron star may rotate during the 22.5 s between SI and S2, and thus four rotation periods may 
be derived for each choice of /^si and ^S2- (The same four times would be derived if instead of 
rotation to S2, ^ = 0, we assume S2, (^ = tt.) 
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Table 1. Cyclotron Scattering Model Parameter Values 



Parameter 


Values 


Bi2 (G) 


(1.55,1.83,2.10,2.40)^ 


^ 


v^'6'3'2'' 


log iVe,21 (cm-2) 


(-2.20, -1.50, -1.10, -0.5, -0.2, 0.1, 0.4, 0.7, 1.0) 


^ 


Bin Size Afi = 0.125 (/i G [0,1]) 


cj) (rad) 


Bin Size A(/> = | 




((/> G [0,7r] if * 7^ or f, 




or ,/. G [0,f ] if ^ = f ) 



^We can shift B12 by it 10% during fits with no loss in model accu- 
racy. 
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Table 2. Frequentist Model Comparison 



1-0 1-1 

Spectrum s^ Polar Cap s^ General Oi^-^MLR ^m Polar Cap s^ General a^^MLR 



SI 


26.7 


(33)* 


23.6 


(31) 


0.21 


37.7 


(33) 


24.4 


(31)* 


1.3x10-3 


S2 


33.6 


(33)* 


31.7 


(31) 


0.39 


33.9 


(33)* 


32.9 


(31) 


0.61 


S1+S2* 


70.8 


(69)* 


59.6 


(66)* 


0.01 


73.7 


(69)* 


59.9 


(65)* 


8.0x10-3 



Note. — The number of degrees of freedom is given in parentheses. 

*Selected model. For both the 1-0 and 1-1 geometries, the data do not conclusively select either 
the polar cap model or general model. 
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Table 3. SI: Bayesian Credible Intervals 



Selected Model: 
Parameter 



la 



1-0 

Polar Cap 

2a 



3cr 



la 



1-1 
General (^ 
2a 



7r\a 
2^ 



3a 



Bi2 (G) 


1 7fi+0-16 1 7fi+0.25 


i-76lS:i 


1.83lgi? 


1 OO+0.28 
^•0-3-0.25 


1 oo-(-0.45 


logiVc,2i (cm-2) 


i.oo+n?'''^ i-oo+^-!i^ 


1.00+°:?° 


1 c;9+0.75 
-i.O^_g2g 


1 c;9+1.58 


-i-52;r7^ 


^ 


U.yyU_QggQ U.yyU„QggQ 


0.990;n^^^ 


n nf^-H0.05h 

'^•uo-o.oe 


0.06+°:Je^ 


0.06+°:^^ 


(f> (rad) 






r, op+0.20i 


n on+o.ee 


o.2o+J:°S 



'^The credible region does not include niarginalization over ^ 

logA''c,2i = 1 is a parameter boundary 
'^The la credible region does not include logiVc^2i = (0.32,0.47) 
Vbest-fit G [0.988,0.992] 

''The lo" credible region does not include fi = (0.09,0.61) 
^The 2a credible region does not include /U = (0.18,0.31) 
^The 3a credible region does not include ^ = (0.187,0.194) 
Vbcst-fit G [0.000,0.125] 
>bcst-fit G [0,0.39 (f )] 
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Table 4. S2: Bayesian Credible Intervals 



1-0 1-1 

Selected Model: Polar Cap General (^ = 0) 

Parameter la 2a 3a la 2a 3a 

B,,iG) 1.96in^ 1-961^:1^ 1-961SJJ l-Sel^e' 1-961^:1^ 1-961^:1 

logiV,^2i (cm-2) -0.22t°-i -0.22+O-i .o.22+?ii -0.22+0-| -0.22+0-^2 -0.22;j-01 

Qoi+0.09a noi+0.49b n Sl+°-6'^ 06+°-^^'= 06+°-32 06+°-^^ 

^ U.Oi_Q 22 U.0i_Q3]^ U.0i_Q3]^ U.UD_QQg U.UD_o o6 U.UD_o o6 



>bcst-fit G [0.250,0.375]. 

'^The 2a credible region does not include fi = (0.56,0.58). 

"^best-fit G [0.000,0.125]. 
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Table 5. S1+S2: Frequentist Model Fits 



P^ s^ Polar Cap 



1-0 



1-1 



s^ General s^ Polar Cap 



s^ General s^^ Paper I^ 



3 


70.8* 


4 


66.5 


5 


64.4 


6 


60.3 



- 


73.7* 


- 


72.2 


69.5 


71.9 


59.6* 


71.3 



71.6 
65.0 
59.9* 



58.0 



* Selected model. 

*The number of free parameters. 

^Using exponentiated Gaussian absorption line profiles (Paper I). 



42 



Table 6. S1+S2: Bayesian Credible Intervals for 1-0 Geometry 







Polar Cap 






General (^ = f )^ 




Parameter 


la 


2a 


3cj 


la 


2a 


3a 


Bi2 (G) 


1 qr;+0-05 


1 qr;+0-ii 
-■-•^^-0.23 


1 qr;+0-15 


1 Q7+0.08 


-I 07+0.14 


-, 07+0.20 

J-^'-0.19 


logiVc,2i (cm-2) 


n on+0.25 
-U.ZZ_Q45 


r, r,r, + 0.55 

-U.^Z_Qg2 


-0.22;?!^ 


r, 00+0.44 


n 99+0.58 
-U.ZZ_Qg2 


-0.22+?:^7 


/^ 


n qi+0.50b,c 


'^•'3J^-0.31 


n q-|+0.69c 








/"SI 








o.i9+°:J« 


o.i9l°:?I 


0.191°:?^ 


MS2 








0.69+HS^ 


0.69i°i^ 


0.69i°i^ 


(/> (rad) 








0.20+°i^^^ 


o.2o:^°|^ 


0.20^°!^ 



'^The credible region does not include marginalization over \I/ 

Vbcst-fit G [0.250,0.375] 

'^The la credible region does not include ^ = (0.49,0.60) 

'^The 2a credible region does not include ^ = (0.94,0.98) 

*^Tlie 3a credible region does not include ^ = (0.97,0.98) 

Vsi,bcst-fit G [0.125,0.250] 

^MS2,bcst-fit G [0.000,0.125] 

Vbest-fit G [0.00,0.39 (f )] 
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Table 7. 81+82: Bayesian Credible Intervals for 1-1 Geometry 







Polar Cap 






General (^ = f )^ 




Parameter 


Icr 


2cr 


3cj 


lo- 


2o- 


3(7 


Bi2 (G) 


1 Qfi+0-06 


i.yD_o.i4 


1 qfi+019 


1 Q/i+0.05 


1.9410:1? 


i.94i°:i 


logiVc,2i (cm-2) 


-U.0Z_Q3g 


n t;o+0.87 
-U.OZ_Q53 


-U.OZ_Qgg 


n 92+0-17 


-0.22+?i8 


-0.22+?:i 


/^ 


n Qf^+0.04b 


'-'•'-"3_o.06 


'-'•'-"3_o.06 








^Si 








n iq+o.o2c 


n iq+0.08 


n iq+0-i6 

U--Ly-0.19 


^S2 








n oi+0.07d 
U-'3J--0.13 


n S1+0-10 


U-3-L-0.31 


(/>Si (rad) 








0.20l°i^'= 


0.20+0|5 


o.2o+J:°^ 


(/>S2 (rad) 








0-59l°:ff 


o.59+°:3'^ 


o.59+°:i 



'^The credible region does not include marginalization over ^ 

Vbest-fit e [0.000,0.125] 

'^si,bcst-fit G [0.125,0.250] 

Vs2,bcst-fit G [0.250,0.375] 

•^^si.best-fit G [0.00,0.39 (f )] 

f0S2,best-fit G [0.39 (f ),0.78 (f )] 
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Table 8. S1+S2: Limits on Neutron Star Rotation Periods 

1-0 (^ = f ) 1-1 (^ = f ) 

lo" 2a 3a la 2a 3a 



irot.l (s^ 


) 2821^^9 


282+^f 


282+^3 


717+1160 
'^'-350 


7,7+10200 


717+00 
'J^'-654 


trot,2 (s^ 


) 65.9;^i 


65.9t\li 


65.9tlli 


55.01^^,^ 


55.0+l'i 


r;c; n+268 


irot,3 (s^ 


) 34.1+5.8 


34.lti^ 


34.lt|^ 


38.0t?i 


SS.O+^s' 


38.0ttl-l 


trot,4 (s^ 


) 24.4+°:^^ 


2AAt\i 


24.4+?:^ 


23.2t'oi 


23.2+2;9 


23.2+J^7^ 
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Table 9. Number of Model Parameters in Fits 





Exponentiated 








Gaussian 


Cyclotron 


Cyclotron 


Spectrum 


Absorption*^ 


Absorption 


Scattering 


SI 


2 


3(1)" 


3 (1-0) 
5 (1-1) 


S2 


2 


5(3)"^ 


3 (1-0) 
3 (1-1) 


S1+S2 


4 


7(3)b 


3 or 6^(1-0) 
3 or 7^(1-1) 



'^See Paper I for description. 

*^Tlie number in parentheses represents the number of 
physically-relevant parameters. 

^The data do not conclusively select either the polar cap 
model or general model. 



1-0 



1-1 



NH 



m 




^^-r> 






Continuum Photon Source 



n 




(a) 

A 



W////////////A 



B 



(b) 




10 20 30 50 70 100 

Energy (keV) 



(a) 

A 



W^MWa 



^mmm^ 



B 



(b) 




10 20 30 50 70 100 

Energy (keV) 



A 



(a) 



n 



B 



(c) 



A 




I I I I I I I I I 1 1 I I I Jfl I I I 1 1 I I I 



— Xi 



1 1 1 1 1 1 1 








o 
o 
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o 



o 

CO 

o 

02 



> 

u 
cu 

5=1 



I 



ifi 




TT 



1_L 



TTT 



J_LL 



O 
O 

in 



o 

O 




O '-I 



I 



ID 



CV2 
I 



CO 
I 



I 



iD 
I 



1- 



A^:^ uo^oiid §OT 



Model Selection 



1 = 1 



M: 



1,1+1 



a 



i,i+2 



a 



i,i+N 




Yes 



Select Simplest 
Alternative Model: 

NewM. 




Yes 



Select M : 



Select M 



N 
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SI 



S2 
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CO 
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oiOOO 
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I ^ ^ L 
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JTjtJ 
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Time (s) 
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60 
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